Free surface instability in a confined suspension jet 
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A jet of non-Brownian particles confined in a thin cell and driven by gravitational force is studied both 
numerically and theoretically. We present a theoretical scheme aimed to describe such a system in the Stokes 
regime. We focus on the dynamics of the interface between the suspension and the pure fluid. Numerical 
simulations solving Newton's equations for all particles show that the jet free surface becomes unstable: the 
fastest growing modes at small sizes coarsen up to the largest structures reaching the jet lateral scale. In the bulk, 
structural waves develop and travel at slightly slower speed than the jet average fall. An analytical model, based 
on hydrodynamic-like equations for the suspension, is derived and predicts the development of the interfacial 
instability. It captures in essence, the collective effects driving the interface destabilization i.e. the long range 
hydrodynamic interactions coupled with the abrupt interface, and no analogous to a surface tension is found. 

PACS numbers: 47.15.Pn, 47.55.Kf, 83.80.Hj 
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I. INTRODUCTION 

Understanding the dynamics of 3D non-Brownian suspen- 
sion in the low Reynolds number regime is a long lasting 
and difficult issue. Long-range hydrodynamics forces cre- 
ate a complex collective particle dynamics 1,2 and up to now, 
no rigorous closed-form formulation of the problem exists at 
moderate densities. For example, difficulties remains such 
as to to explain particle sedimentation, dispersion and mix- 
ing in a finite size container (see 3 and references inside). Re- 
cently, the miniaturization of hydrodynamics devices^ nec- 
essary to develop microfluidic apparatus performing at low 
Reynolds number, mixing or separation, has brought to the 
front the importance of boundary confinement of suspensions. 
In fact, numerous devices bear the form of narrow channels or 
Hele-Shaw cells with one direction reaching only few parti- 
cle sizes&. Note that experimentally quasi-2D non-Brownian 
suspensions have also been studied in the case where the cell 
gap was almost the particle size 7 - 8 - 9 . In this case, the existence 
of anti-drag correlations due to fluid recirculation around each 
grains were identified 9 . 

To get a better description of mixing, an important basic 
issue is to capture the evolution of an interface between a sus- 
pension and a particle-free fluid (the pure fluid). This problem 
is also related to many studies done to understand miscible in- 
terfaces dynamics, either from two different fluids 10 or a fluid 
in contact with a suspension falling under gravityAA*A2ii^. 

Here, we describe the evolution of a jet of suspended parti- 
cles in a thin cell driven by gravitational forces. The suspen- 
sion is non Brownian and the hydrodynamic forces between 
particles are obtained in the Stokes regime. Numerical simu- 
lations that solve the Newton equations for all particles follow 
the evolution of the jet free surfaces. We present an analytical 
model, based on hydrodynamic-like equations for the suspen- 
sion which is able to predict the development of the instability. 
We present and extend a scheme that was developed recently, 
in the context of a 2D suspension cluster falling downil and 
also presented briefly in ref«I£ in the context of a jet. The rel- 



ative simplicity of the model, focusing on the actual effect of 
long range hydrodynamic forces, allows to discuss the domi- 
nant physical features of the interface dynamics. 



II. MODEL 

We consider a system of N solid particles that move through 
an incompressible Newtonian fluid of viscosity 77. The fluid is 
confined between two parallel plates separated at a distance 
2d in the z direction, being infinite in the other two directions. 
To simplify the computation of the hydrodynamic interaction 
forces, and to focus on the effect of the long-range nature of 
them, we consider cylindrical particles of height L (slightly 
smaller than 2d) and radius cr that have planar motion only. 
Particles are thin, i.e. 2d <K cr and their mass is m. 

The hydrodynamic forces between the cylinders in the con- 
fined geometry have been computed in RefJi in the Stokes 
regime for a dilute suspension. First, the force over the z'-th 
particle has a drag component given by -^u/, where m, is the 
in-plane velocity of particle i and t\ = md/ncr 2 ?] is the relax- 
ation time of a single particle. Also, the presence of the other 
particles induce hydrodynamic forces over each other. When 
particles are far apart, the force on z due to the presence of k is 
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(1) 



where Rki = Ri - Rk is the relative distance between particles, 
being R, the position of the center of mass of the z'-th particle 
and the tensor K is given by 



K(R) = (cr/R) 2 (I - 2RR/R 2 ) 



(2) 



This force depends both on the direction of the velocity zt* 
and on the relative distance R^. When z4 is parallel to R%, 
the interaction force on i is parallel to z4 (drag), and if Uk is 
perpendicular to the force turns out to be in an opposite 
direction to z4 (antidrag). 



2 



On the other hand, when particles are close to each other, 
lubrication forces appear producing the net force on particle i 
due to the presence of k— 

Pf k = 2ndr 1 ^-^ I {R 2 ik l-R ik R ik )il ik 

V e R ik 

3 d 2 1 Uik- Rik g 

+ -2ndii—r -r — R ik (3) 

2 cr £ e R z . 

ik 

where e = Rtk/cr — 1 is the gap between the particles. The 
force depends on the relative velocities between the pair and 
it respects the action-reaction principle. Therefore it does not 
change the total momentum of the pair, but reduces the rela- 
tive velocity. 

The cutoff for using either expression Q or Q is rather 
arbitrary. We adopt the convention that when the pair is closer 
than /?i u br, the lubrication force is used, when the distance 
is larger than Rf ar , the far force Q is used, and in between a 
linear interpolation between the two is computed. The result 
is the interaction force PL. 

There is also an extra drag force produced by the flow be- 
tween the particles and the plates 14 . This force can be added 
to the other drag force, simply modifying the prefactor by a 
value given by the particular experimental setup. In order to 
simplify the analysis, in what follows we will disregard the 
presence of this term but in case of being relevant to some 
experimental configuration, it can be trivially included. 

In summary, the dynamical equations for the suspended 
particles are 

dui m v-i -», 

m— = -— Ui + 2_^F ik + mg (4) 

that constitute a complex system of equations that will be 
solved numerically. 

Although the tensor diverges at close distances, it has 
the remarkable property that, its integral over the direction of 
R vanishes. This implies that if a particle is surrounded by an 
homogeneous distribution of suspended particles, all moving 
with equal velocities, the far-field contribution to the hydro- 
dynamic force vanishes. 

The long range nature of the forces forces between particles 
makes it computationally expensive to perform direct numer- 
ical simulations where all pairs of forces are computed. How- 
ever, the expressions Q and (0 for the far-field contribution 
to the hydrodynamic forces, suggests that a mean field approx- 
imation can be done. In fact, if particles in a region move with 
similar velocities they exert similar forces to a target particle; 
it is natural then to group them altogether and compute the to- 
tal force made by the group on any target particle. We define 
the coarse grained current J on a cell c as 

/c = ay2"-' (5) 

fee 

where AS is the area of each cell. 



Then, if the system is nearly homogeneous in each cell, the 
far field contribution to the force (0 can be approximated by 

^i = ~ X AS K(Ri - Rc)Jc (6) 

c 

where the sum runs on all cells, R c is the position of the cen- 
ter of each cell. This expression allows to study systems of 
many particles at reasonable cost. Note that the vanishing in- 
tegral of K(R) over the angles implies that (|6) can be applied 
down to neighboring cells without producing divergences, but 
the mean field approximation can be inaccurate for near cells. 
Local inhomogeneities are not completely taken into account 
by the mean field and unrealistic vanishing forces can be ob- 
tained. To solve this problem, the mean field approximation of 
the far field force is applied only to cells that are separated by 
a distance larger than R m f, otherwise the direct summation on 
the pair of particles is performed. Finally, the values chosen 
for performing the simulations are R\ u b r = 1.3<x, R^ = 2.0<x, 
R m f = 5.66<x. This numerical method was shown to give accu- 
rate results in the study of the spreading of a falling cluster 14 . 



III. SIMULATIONS OF A FALLING JET 

The long range forces Q vanish when a particle is sur- 
rounded by and homogeneous medium, but in presence of 
inhomogeneities the force is finite. The effect of inhomo- 
geneities is clearly seen when there a separation line between 
a region with suspended particles and a region without them; 
it was shown in 14 , that the force is enhanced if the separation 
line is curved, being possible to produce instabilities. As our 
description is for the dynamics of the suspended particles, and 
not for the surrounding fluid, it is sensible to call free surface 
the separation line between a region with suspended particles 
and the pure fluid. 

In order to study the effect of the long range forces on the 
free surfaces and, in particular, if they can lead to unstable 
motion, we consider a falling jet of particles immersed in a 
fluid. 

A jet of falling particles is studied numerically. Initially, 
a jet of N - 24000 particles is placed randomly at rest in a 
rectangle of width L x = 90cr and height L v = 600<x. They fall 
down due to the action of a gravitational field pointing in the 
positive y direction. To mimic an infinity jet the vertical direc- 
tion is periodic and the force computation uses the minimum 
image convention 16 . Units are chosen such that the particle 
diameter cr, particle mass m, and the relaxation time t\ are set 
to one. The gravitational force is mg =1.0 and therefore the 
limiting velocity for a single particle V™ = gT\ = 1 .0 (see Eq. 

Fig- HI shows three successive snapshots of the jet. They 
show that the free surfaces become unstable, showing the ap- 
pearance of waves. At the beginning the surface waves are 
characterized by short wavelengths but later as the amplitude 
of the perturbations growth, the characteristic wavelengths 
grow also. A coarsening process is developed leading to larger 
structures. Once the size of this structures is comparable to the 
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jet width, interactions between the two surfaces are observed 
and in-phase surface oscillations are obtained. The surface 
waves are also accompanied by modulations of the particle 
density. 
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FIG. 1: Numerical simulation of a suspension of N = 24000 parti- 
cles falling by the action of gravitational field. From the left to the 
right: t = 200, t = 2000, and t = 4000. Units are described in the 
text. 



Fig. |2 shows the average vertical velocity of the falling 
jet as a function of time. It is seen that the jet rapidly gets 
an asymptotic velocity that is larger than the one for a single 
isolated particle. At much larger times, as the jet develops 
large structures and the particles separate, the falling velocity 
decreases and approaches V™ = 1 .0. 




FIG. 2: Instantaneous average of the vertical velocity as a function of 
time (continuous line). The dashed line is the asymptotic theoretical 
value of the jet if it does not deform, V™ = 1.165. Inset: evolution 
for longer times. 



To describe in more detail the coarsening process, the x- 
integrated density N x (y, f) is studied, computed as the number 
of particles in a horizontal bin of height <j centered at y. In 



Fig. |3]a spatio-temporal plot N x (y, t) is presented in the refer- 
ence frame falling with the same instantaneous velocity of the 
jet. It is seen that structures grow in time showing the coars- 
ening process. Also, the structures propagate in the negative y 
direction; as the plot is presented in the reference frame of the 
falling jet it indicates that the perturbations fall at a slightly 
smaller velocity than the jet. 




FIG. 3: Spatiotemporal diagram of the evolution of the x-integrated 
density N x (y, t) in the reference frame of the falling jet. Time is on 
the vertical axis and increasing upward. The horizontal axis is the y 
coordinate, with periodic boundary conditions, and the gravitational 
acceleration points to the right. The gray scale is proportional to den- 
sity, with lighter regions representing denser regions of the system, 
the minimum and maximum represented densities are Nf m = 11 and 
yymax _ ^ respectively, and the average is N^ s = N/L y = 40. The 
simulation parameters are the same as in Fig.Qand the total simula- 
tion time is t = 6152. 



IV. GLOBAL MODEL 

The behavior observed in the numerical simulations suggest 
an instability like the observed when two immiscible fluids 
in contact move with a relative velocity (Kelvin-Helmholtz 
instability )ii 

To describe the appearance of the surface instability, we 
build a global model, similar to hydrodynamic equations. We 
consider the particle number density n(R) and the particle 
mean velocity V(R). The average force density over the sus- 
pension, produced by the drag force, the far force contribution 
Q, and the gravity acceleration is 

i^ i -, r -> -> -> ^ 

F(R) = J(R) n(R) dR'K(R-R')J(R') + n(R)mg 

T\ 8ti J 

(7) 

with J = mnV is the mass current density. 

In a Euler-like global model as the one proposed, the lubri- 
cation force contribution can be neglected because its effect 
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is to reduce velocity fluctuations, but it does not modify the 
mean velocity as it only affects the relative velocity. Even- 
tually, its effect would be to produce an additional effective 
viscosity as obtained in kinetic theoryi™. 

Therefore a global equation of motion for the suspension 
can be written as 



dt v ' 



nm 



(8V -i 

_ + y . w 



= P 



(8) 
(9) 



Before analyzing the development of the instabilities, let's 
consider the evolution of a unperturbed homogeneous jet. For 
a thin jet of width L x and height L v , with periodic boundary 
conditions in y, if the particle current is assumed to be homo- 
geneous J(R, t) = Jj et {t)y, the integral term in simplifies 
(see Eq. iA2i ) and we have that 

J dR'K(R - R')Jj et 

= 2cr 2 [arctan(L v /L. t ) - arctan(L r /L v )] /j et }> (10) 

independent of Therefore Eqs. and admit a solution 
of an homogeneous falling jet with velocity 



V jet (t) = gT }et (l-e- t/r >«) 



where the jet relaxation time is 



'jet 



1 - «ocr 2 [arctan(L v /L t ) - arctan(L x /L v )] /4 



(11) 



(12) 



and no is the jet number density. 

Note that the asymptotic velocity as the relaxation time is 
different from the one of a single falling particle V™ = gT\ . 
This effect is due to the long range hydrodynamic interactions. 
A peculiar feature is that the modification of the falling veloc- 
ity depends only on the system shape, but not on its size. A 
similar phenomena was observed in the case of a falling circu- 
lar clustesil. In fact, for a thin jet (L x < L y ) most of the relative 
distances between particles are parallel to the falling velocity, 
giving rise to a positive hydrodynamic forces between parti- 
cles. That is particles exerts a positive drag, in the direction 
of the motion, increasing the falling velocity. On the other 
hand, if the jet were wide (L x > L v ) most of the relative dis- 
tances between particles would be perpendicular to the falling 
velocity and the dominant contribution of the hydrodynamic 
forces would be the antidrag, reducing the jet falling velocity. 
In the case of the jet we are considering the theoretical values 
for the falling jet velocity and relaxation time are V™ t = 1 . 165 
and Tjet = 1.165, that are larger than the values for a single 
particle. In Fig. [5]it is seen that the jet rapidly gets this ve- 
locity and afterward, when it starts deforming, the velocity 
slowly decreases approaching V™ = 1 .0. This long time value 
is consistent with Eqs. dl It and dl 2t . when the jet density de- 
creases. An exponential fit to the early time evolution of the 
velocity of the jet gives the same relaxation time as predicted. 



n 
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FIG. 4: Left: Sketch of the semi - infinity suspension. Q.„ is the un- 
perturbed domain, fli is the perturbed domain, £ is the perturbation 
in the x direction. 



In presence of a free surface, an additional equation must 
be added to describe the evolution of the free surface posi- 
tion £(j, f) (see Fig. @}. The mean suspension velocity at the 
surface must be equal to the surface velocityii. 



t=f oy 



JC=f 



(13) 



V. LINEAR STABILITY ANALYSIS 



The stability of the free surface of the jet is analyzed using 
the global Eqs. and dl 3I >. To simplify the problem and 

assuming that at the beginning the two free surfaces do not 
interact strongly, we will consider the simpler case of a single 
free surface, limiting a semi-infinite homogeneous suspension 
in the x < region, of density «o- In this geometry the jet 
falling velocity is Vo = gT\y/(l + jSo/4), where fio = 7nx 2 rco/4 
is the area fraction of the suspension (see Eq. dA3» . Note that 
in this geometry Vo < V™. 

We consider small perturbations of the free surface. Intro- 
ducing a formal parameter e <K 1, the fields are written as 
V(R) = Vo + eVdR), n = n + en x (R), and <f = e£i(y). The 
domain O over which the integration in Q must be done can 
be split into the initial domain Qo : x < 0, plus a perturbed 
domain Oi : < x < g(y, t) (see Fig. [4}. Keeping terms up to 
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linear order in e, the equations read 

^7 + « () V • ft = 
at 

dft 1 ^ ^ 1 B, 

— 1 + y • Wi = ViCR) - — -^K(y)J 

at t\ 2ti 4 

-^-tio f dR'K(R - R')fi(R') 

i r°° 

-—no d/K(/?-/y)/ £Cy') 
2ti 

3f I d£ I 

-£ + V M ,T = V A n 

df ' I .r=0 cry b=0 



(14) 



(15) 
(16) 



The deformed surface has a long range effect of the fluid mo- 
tion, as reflected by the last term in Eq. [2] 

This integro-differential system of equations does not have 
an evident analytic solution. However, an estimation of the 
instability modes can be obtained as follows. A Fourier ex- 
pansion of the fields Ai(R, t) = A\{k, f)e is performed. As 
the Fourier basis is not a solution of the system (particularly 
because of the integral terms in (I15» . we obtain an estimation 
by evaluating the equations at x — 0. Defining dimensionless 

variables p = hjrio, U = ft/(gri), if = f/(grf), q = (gT 2 )t 
and s = t /ti , the equations read 



dp _ -> . 

— = -iq- U - iq ■ U p (17) 
us 

— =-&- f} a )C - Qi0+ ptf ) - *h&o<P (18) 

OS 



dtp _, -* 

— = U x - i(q ■ Uo)<p 
os 

with the following defined tensors (see Eq. (|A5jl) 



i 1 



Q ' " 2ll0 -lj + ^ 

Q2 = -\q y \\ • j 



(19) 

(20) 
(21) 



We recall the the dimensionless jet falling velocity is Uq = 

57(1 +yS /4) 

The tensors Qi and Q2 capture the effect of the long range 
forces produced by the perturbations, when they are integrated 
in the semi-infinite volume. Q t describes the interaction of the 
perturbations with the bulk and Q2 the interaction with the free 
surface. Remarkably, Q2 is proportional to q and not to q 2 , 
therefore the effect of a curved surface cannot be described in 
terms of a surface tension. Note also that Qi is independent of 
the magnitude of q and depends only on its direction. This fact 
has an important consequence because in the limit of pertur- 
bations of small wave vectors the effect of the hydrodynamic 
interactions does not vanish. In fact in the limit q — > the 
linear system of equations ( I19> has a Jordan-block structure 
and admits a solution of the form (at first order in /3o) 



p = 2AU '(- cos (p + i\ sin< 

U x = A>8ol sin 0| 

U y = -Afiocoscf) 

ip = Afio\ sin <f>\t 



(22) 
(23) 
(24) 
(25) 



where A is an arbitrary coefficient given by the initial condi- 
tion and 4> is the angle between the wave vector q and the x 
direction. The solution shows a linear, instead of exponential, 
increase in time of the surface oscillations. 

For larger values of q the system of equations looses its 
Jordan-block structure and solutions in the form exp(As) are 
looked for. The eigenvalues A for q <k 1 , ySo <k 1 , and < (p < 
7T (sin d> > 0) are 



x x = 


-iU q(l +/3o/4) sin0 


V3 

+ — p {) U qs\n(, 


b + 0{q 2 t)26) 


A 2 = 


-iU q(l +/?o/4) sin0 


V3 

- —foUoq sin (, 


b + O{q 2 021) 


X 3 = 


-1 — —(cos(p + 1 sin. 


P) + 0(g) 


(28) 


A 4 = 


1 #>/ ± ■ ■ 
-1 + — (cos0 + 1 sin. 


t) + O(q) 


(29) 



and similar results for n < (p < 2n (sin <p < 0). 

Two of the eigenvalues, A3 and A4 have negative real parts 
for small q and therefore correspond to damped motion. How- 
ever, the real part of A\ (A2 for n < <p < 2n) is positive. There- 
fore an instability is predicted. The analysis shows that the 
system becomes unstable for any strength of the gravitational 
force and, coming back to the original units, the instability 
rate is directly proportional to Vok. In the limit k «c 1 the 
real parts of A[ and Xi are proportional to the wavevector. It is 
reasonable to expect that a more detailed model, that includes 
terms proportional to gradients of V or the viscous effect due 
to the lubrication forces, will produce that for high enough 
wave vectors, the real part of the eigenvalues become nega- 
tive again. Hence, it is expected that the system is unstable 
for a range of wave vectors going from zero to a finite value, 
showing the coarsening process observed in the simulation 19 . 

The eigenvalues A[ and A2 have also an imaginary part. 
The dominant one with respect to k and /?o is A\j = 
-iqUo(l + y6o/4)sin0. It corresponds to a wave in the form 
exp(iq sin (p(y - Uq(1 + /3o/4)t) + iq cos <px), that propagates 
in the +y direction, with phase velocity C/o(l + ySo/4) = 1. 
That is, the structures propagate at the same velocity than a 
single falling particle. This phase velocity is larger than the 
jet falling velocity by a factor 0(J3q). As a consequence, the 
unstable structures should move faster than the jet. In the sim- 
ulations, we actually find the opposite: the structures move 
slower than the jet. This change of character can be due to the 
different geometries as they are considered, as it also happens 
to the jet falling velocity that in one case is faster than V™, 
while in the other case was found slower. 

The linear stability analysis for a jet with two free surfaces 
as in the simulations is much more involved. One extra equa- 
tion must be added and the long range interactions in this ge- 
ometry give more complex expressions, as can be seen in the 
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case of the falling velocity. Nevertheless, we expect that the 
analysis presented here, showing that the surface is uncondi- 
tionally unstable and the growth rates are proportional to the 
wave vectors, is preserved. 

We have studied numerically the growth rate for different 
wavevectors in the jet. For each wavevector we have per- 
formed new simulations, with an initial condition similar to 
the one described previously for the jet except that the x co- 
ordinates are modulated in such a way to create two waves 
of vectors k y = 2im y jL y in the surface positions. In prac- 
tice we map the initial x coordinates of the original rectangu- 
lar jet as [0, L x ] — > [-A cos(k y y), L x + A cos(k y y)]. In each 
simulation we compute the evolution of the Fourier ampli- 
tudes = AT 1 ! £,.(>,■ - L A ./2) 2 e 2m '">' /L > | = AT 1 1 J p(x,y)(x- 
L x /2) 2 e 2mny t L >'\, with the same wavenumber as the perturba- 
tion. These Fourier amplitudes measure, at linear order, the 
amplitude of g(ky). Fig. |5] shows the results of the Fourier 
amplitudes for the simulations done with n y = 1,2,3,4. For 
n y = 1 a linear increase of is observed while for n y > 2 
exponential increase of the *P's are obtained, until nonlin- 
ear interaction between modes appear. The growth rates are 
A(n y = 2) = 1.21 x 1(T 3 , A(n y = 3) = 1.82 x 1(T 3 , and 
A(n y = 4) = 2.19 x 1CT 3 , that are roughly in the ratio 2:3:4, 
confirming the linear proportionality with k. However, no di- 
rect comparison with the prediction can be made because only 
k y was fixed but not k x , therefore the angle <p takes all possi- 
ble values. Larger values of n v give very noisy results to be 
cor 

0.01 1 , 1 , 1 y , 

0.008- 

0.006- '' /' 




°0 500 1000 1500 2000 

t 

FIG. 5: Evolution of the Fourier amplitudes *?„ computed with the 
same wavenumber n as the perturbation made on the surface, n = 1 
continuous line, n = 2 dotted line, n = 3 dashed line, and n = 4 
dash-dotted lined. 



VI. CONCLUSIONS 

In conclusion, we investigated the dynamical evolution, in 
the Stokes regime, of a jet of falling particles confined in an 
Hele-Shaw cell. The free surfaces of the jet is shown to be- 
come unstable and a coarsening process takes place from the 
smallest sizes up to the largest structures. 

A continuous Euler-like hydrodynamic model for the sus- 
pension is presented which describes the main features of the 
jet flow. In particular, the model predicts a geometry depen- 
dent falling velocity for the jet (thin jets fall slower than thick 
jets). 

The theoretical analysis shows that the free surface is unsta- 
ble to any gravitational force and any wave vector. The growth 
rates are proportional to the wave vectors, result that is consis- 
tent with the coarsening process observed in the simulations. 
Finally, the model predicts that the structures move with a ve- 
locity slightly different than the jet's velocity, as observed in 
the simulations but with a different sign prediction. In this 
problem we find that the interface dynamics is quite different 
than the usual immiscible fluid case where the surface tension 
play a leading role such as to select and stabilize structures. 
For the present jet,the particle hydrodynamic interactions cre- 
ate a surface instability essentially due to long range effects 
caused by flow recirculation around a particle. A crucial point 
is that this effect is scale free. We observe at the end of the 
simulations a partial mixing due to flow structuration as well 
as hydrodynamic dispersion causing a diffuse interface. 

Finally, we open a question on the relation between this 
model and the dynamics of actual non-Brownian suspensions 
in a confined cell. Though we have here a situation a little bit 
artificial as the problem is solved in the case of particle sizes 
(diameter) larger than the cell gap, we believe that all the fea- 
tures dealing with scales larger than the gap ( essentially the 
long range interactions) will survive at the qualitative level. 
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where R = xx + yy is any vector in the interior of the integra- 
tion domain and we have used periodic boundary conditions 
in y. When L y » L y the integral depends slightly on x and its 
horizontal average is in that limit 



h » 2cr 2 [arctan(L v /L x ) - arctan(L t /L v )] | * Jj* j (A2) 



(b) Homogenenous integral for the semiinfinite system 



X0 poo TTCp' l 1 \ 

dx' J dy'K(M-M') = -—[ L (A3) 

with 3 — xx + yy any vector with x < and we have used 
periodic boundary conditions in y. 

(c) Fourier transform for the semiinfinite system 



X0 poo 
dx I dy'K(R-R') 
CO %J — oo 



- ncr 



\ky 



1 

0-1) ik x + \k y \ \-i -1 



JkR' 



1 -i 



(A4) 



with R — xx + yy any vector with x < and we have used 
periodic boundary conditions in y. 

(d) Fourier transform for the interface line 



APPENDIX A: KERNEL INTEGRALS 

Some integrals of the kernel K(R) used in the 
manuscript are: 

(a) Homogenenous integral for the jet 



h = f dx f dy' K(R - R') 

J-L, J-Ly 



2cr z 



arctan 



Ly X 



(L x + x 
arctan 



1 

-1 



(Al) 



h 



-f 



dy K.(xx + yy)e 
1 



iky 



(A5) 



The Dirac delta appears when K is evaluated up to the ori- 
gin, that physically cannot happen because the interaction is 
replaced by the lubrications forces. Therefore, this term will 
be eliminated in the evaluation of the kernel integrals. 



